Langesim basic usage
[1]:
from langesim import Simulator
Define the stifness k(t) and center c(t) of the harmonic potential
[2]:
def k(t):
"""Stifness of the potential.
Varies linearly from ki=0.5 to kf=1.0 in tf=5.0 time
units, then remain constant at kf
"""
tf = 5.0
ki = 0.5
kf = 1.0
if t <= 0.0:
return ki
if t < tf:
return ki + (kf - ki) * t / tf
if t >= tf:
return kf
def c(t):
"""Center of the harmonic potential.
Varies linearly from ci = -1.0 to cf = 1.0 in tf time units, then remains constant.
"""
tf = 5.0
ci = -1.0
cf = 1.0
if t <= 0.0:
return ci
if t < tf:
return ci + (cf - ci) * t / tf
if t >= tf:
return cf
Create the simulator
[3]:
simulator = Simulator(
tot_sims = 100_000,
dt = 0.001,
tot_steps = 10_000,
snapshot_step = 100,
k = k,
center = c
)
The simulator has not performed any simulations yet:
[4]:
print(simulator.simulations_performed)
0
Run a simulation with the same parameters given in the construction of the simulator
[5]:
simulator.run()
One simulation has been performed
[6]:
print(simulator.simulations_performed)
1
The simulation can be accessed at simulator.simulation[0]
[7]:
simulator.simulation[0].results
[7]:
{'times': array([ 0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1. ,
1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2. , 2.1,
2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9, 3. , 3.1, 3.2,
3.3, 3.4, 3.5, 3.6, 3.7, 3.8, 3.9, 4. , 4.1, 4.2, 4.3,
4.4, 4.5, 4.6, 4.7, 4.8, 4.9, 5. , 5.1, 5.2, 5.3, 5.4,
5.5, 5.6, 5.7, 5.8, 5.9, 6. , 6.1, 6.2, 6.3, 6.4, 6.5,
6.6, 6.7, 6.8, 6.9, 7. , 7.1, 7.2, 7.3, 7.4, 7.5, 7.6,
7.7, 7.8, 7.9, 8. , 8.1, 8.2, 8.3, 8.4, 8.5, 8.6, 8.7,
8.8, 8.9, 9. , 9.1, 9.2, 9.3, 9.4, 9.5, 9.6, 9.7, 9.8,
9.9, 10. ]),
'x': array([[-3.65804544, -3.13428626, -2.51367711, ..., 1.31865933,
1.42542587, 1.60251312],
[-2.03745209, -2.37827146, -1.84574698, ..., 1.29004746,
1.5636284 , 1.52115507],
[-3.46189917, -3.2660894 , -3.56129824, ..., 0.48226885,
0.4256303 , 0.95760999],
...,
[-1.17028482, -1.39109569, -1.17040899, ..., 0.83817809,
1.09375208, 0.63111152],
[-1.00796136, -0.40103226, -0.47672899, ..., 2.30937581,
2.34016721, 2.18571268],
[-1.50659137, -2.28013269, -1.69045596, ..., -0.22986617,
0.16568942, 0.78374506]]),
'power': array([[ 0. , 0.67980267, 0.45836984, ..., 0. ,
0. , 0. ],
[ 0. , 0.38980455, 0.23532712, ..., 0. ,
0. , 0. ],
[ 0. , 0.73621162, 0.89806561, ..., 0. ,
0. , 0. ],
...,
[ 0. , 0.09717766, 0.05516869, ..., 0. ,
0. , 0. ],
[ 0. , -0.09842561, -0.08239977, ..., 0. ,
0. , 0. ],
[ 0. , 0.35635099, 0.18986255, ..., 0. ,
0. , 0. ]]),
'work': array([[ 0.00000000e+00, 7.77702863e-02, 1.29297753e-01, ...,
3.50858817e+00, 3.50858817e+00, 3.50858817e+00],
[ 0.00000000e+00, 3.00013965e-02, 5.67694600e-02, ...,
2.11392028e+00, 2.11392028e+00, 2.11392028e+00],
[ 0.00000000e+00, 7.50189433e-02, 1.61997221e-01, ...,
9.81078639e-01, 9.81078639e-01, 9.81078639e-01],
...,
[ 0.00000000e+00, 6.39956596e-03, 1.11596362e-02, ...,
1.25889066e+00, 1.25889066e+00, 1.25889066e+00],
[ 0.00000000e+00, -1.54242454e-03, -9.45528345e-03, ...,
1.81095009e+00, 1.81095009e+00, 1.81095009e+00],
[ 0.00000000e+00, 2.68537198e-02, 5.79637727e-02, ...,
1.58282763e+00, 1.58282763e+00, 1.58282763e+00]]),
'heat': array([[ 0.00000000e+00, -6.38553886e-01, -1.23524940e+00, ...,
-5.22411768e+00, -5.18439598e+00, -5.09337853e+00],
[ 0.00000000e+00, 2.13852846e-01, -1.03024228e-01, ...,
-2.34093322e+00, -2.22415851e+00, -2.24719569e+00],
[ 0.00000000e+00, -2.34153494e-01, 1.36644564e-01, ...,
-2.36229274e+00, -2.33136524e+00, -2.49541706e+00],
...,
[ 0.00000000e+00, 3.37412952e-02, -2.10565429e-03, ...,
-1.25304672e+00, -1.26174516e+00, -1.19810053e+00],
[ 0.00000000e+00, 8.12000368e-02, 6.05266255e-02, ...,
-9.53733439e-01, -9.12941868e-01, -1.10800866e+00],
[ 0.00000000e+00, 3.53388909e-01, 3.22141444e-02, ...,
-8.90700933e-01, -1.29894927e+00, -1.62360323e+00]]),
'delta_U': array([[ 0. , -0.5607836 , -1.10595164, ..., -1.71552951,
-1.67580781, -1.58479036],
[ 0. , 0.24385424, -0.04625477, ..., -0.22701294,
-0.11023822, -0.1332754 ],
[ 0. , -0.15913455, 0.29864178, ..., -1.38121411,
-1.3502866 , -1.51433842],
...,
[ 0. , 0.04014086, 0.00905398, ..., 0.00584394,
-0.0028545 , 0.06079013],
[ 0. , 0.07965761, 0.05107134, ..., 0.85721666,
0.89800823, 0.70294143],
[ 0. , 0.38024263, 0.09017792, ..., 0.6921267 ,
0.28387837, -0.0407756 ]]),
'energy': array([[1.76630139e+00, 1.20551779e+00, 6.60349751e-01, ...,
5.07718845e-02, 9.04935850e-02, 1.81511030e-01],
[2.69076708e-01, 5.12930950e-01, 2.22821940e-01, ...,
4.20637659e-02, 1.58838484e-01, 1.35801304e-01],
[1.51523688e+00, 1.35610233e+00, 1.81387866e+00, ...,
1.34022771e-01, 1.64950275e-01, 8.98456672e-04],
...,
[7.24922968e-03, 4.73900909e-02, 1.63032116e-02, ...,
1.30931652e-02, 4.39472623e-03, 6.80393552e-02],
[1.58458088e-05, 7.96734581e-02, 5.10871878e-02, ...,
8.57232502e-01, 8.98024073e-01, 7.02957280e-01],
[6.41587029e-02, 4.44401332e-01, 1.54336620e-01, ...,
7.56285402e-01, 3.48037069e-01, 2.33830998e-02]])}
Let’s perform another simulation with a different set of parameters.
[8]:
simulator.run(
tot_sims = 100_00,
dt = 0.01,
tot_steps = 7_000,
snapshot_step = 500
)
Now 2 simulations have been performed
[9]:
print(simulator.simulations_performed)
2
The second simulation can be accessed at simulator.simulation[1]
[10]:
simulator.simulation[1].results
[10]:
{'times': array([ 0., 5., 10., 15., 20., 25., 30., 35., 40., 45., 50., 55., 60.,
65., 70.]),
'x': array([[ 0.6518781 , 1.65103553, -0.58825238, ..., 0.32742353,
0.23934235, 3.05964769],
[-2.28417022, -0.27663891, 0.34159657, ..., 0.17605995,
1.59999249, -0.1653596 ],
[-2.09288736, -0.3504792 , 0.06292072, ..., -0.14239458,
-0.32201354, 1.17601181],
...,
[-1.63223156, 1.50262862, 0.51653385, ..., 1.07197109,
1.48655274, 0.40872803],
[-2.88004875, -0.5673787 , 2.66704629, ..., -0.04229329,
2.50997465, -0.42593398],
[-2.66714664, 0.43006473, 0.95055973, ..., -1.15124269,
-0.54634716, 0.18006437]]),
'power': array([[ 0. , -0.23976064, 0. , ..., 0. ,
0. , 0. ],
[ 0. , 0.59083605, 0. , ..., 0. ,
0. , 0. ],
[ 0. , 0.63004199, 0. , ..., 0. ,
0. , 0. ],
...,
[ 0. , -0.18901782, 0. , ..., 0. ,
0. , 0. ],
[ 0. , 0.74835913, 0. , ..., 0. ,
0. , 0. ],
[ 0. , 0.24318825, 0. , ..., 0. ,
0. , 0. ]]),
'work': array([[ 0. , -0.19768468, -0.19768468, ..., -0.19768468,
-0.19768468, -0.19768468],
[ 0. , 2.17710764, 2.17710764, ..., 2.17710764,
2.17710764, 2.17710764],
[ 0. , 3.41178549, 3.41178549, ..., 3.41178549,
3.41178549, 3.41178549],
...,
[ 0. , -0.491979 , -0.491979 , ..., -0.491979 ,
-0.491979 , -0.491979 ],
[ 0. , 1.62537502, 1.62537502, ..., 1.62537502,
1.62537502, 1.62537502],
[ 0. , 0.68294696, 0.68294696, ..., 0.68294696,
0.68294696, 0.68294696]]),
'heat': array([[ 0. , -0.272567 , 0.77678218, ..., -0.25831108,
-0.19519061, 1.63658367],
[ 0. , -1.77447748, -2.37263339, ..., -2.24994233,
-2.40938543, -1.91034943],
[ 0. , -2.79848916, -3.2713274 , ..., -3.0578535 ,
-2.83652629, -3.69489611],
...,
[ 0. , 0.51836758, 0.50891958, ..., 0.39463973,
0.5104166 , 0.56685109],
[ 0. , -1.28068284, -1.11949918, ..., -1.9658332 ,
-1.36900912, -1.49237699],
[ 0. , -1.21537833, -1.37656927, ..., 0.93613111,
-0.18219667, -1.04164421]]),
'delta_U': array([[ 0. , -0.47025168, 0.5790975 , ..., -0.45599576,
-0.39287528, 1.43889899],
[ 0. , 0.40263016, -0.19552576, ..., -0.07283469,
-0.23227779, 0.2667582 ],
[ 0. , 0.61329634, 0.14045809, ..., 0.35393199,
0.5752592 , -0.28311062],
...,
[ 0. , 0.02638858, 0.01694057, ..., -0.09733927,
0.0184376 , 0.07487209],
[ 0. , 0.34469218, 0.50587584, ..., -0.34045818,
0.2563659 , 0.13299803],
[ 0. , -0.53243137, -0.69362231, ..., 1.61907807,
0.50075029, -0.35869726]]),
'energy': array([[6.82175310e-01, 2.11923634e-01, 1.26127281e+00, ...,
2.26179551e-01, 2.89300027e-01, 2.12107430e+00],
[4.12273291e-01, 8.14903449e-01, 2.16747536e-01, ...,
3.39438601e-01, 1.79995496e-01, 6.79031493e-01],
[2.98600695e-01, 9.11897032e-01, 4.39058786e-01, ...,
6.52532687e-01, 8.73859895e-01, 1.54900784e-02],
...,
[9.99291876e-02, 1.26317764e-01, 1.16869761e-01, ...,
2.58991879e-03, 1.18366786e-01, 1.74801274e-01],
[8.83645824e-01, 1.22833800e+00, 1.38952166e+00, ...,
5.43187646e-01, 1.14001172e+00, 1.01664385e+00],
[6.94844478e-01, 1.62413108e-01, 1.22216999e-03, ...,
2.31392255e+00, 1.19559476e+00, 3.36147221e-01]])}
Let’s analyse the results of the first simulation
[11]:
sim = simulator.simulation[0]
Plot trajectories number 1, 2, 5 and 10
[12]:
sim.plot_sim('x', [1, 2, 5, 10])
One can plot the following quantities:
[13]:
sim.result_labels
[13]:
['x', 'power', 'work', 'heat', 'delta_U', 'energy']
Plot the work done of trajectories 0 to 9
[14]:
sim.plot_sim('work', range(10))
Plot the average position, the variance of the position, the average work and heat
[15]:
sim.plot_average('x')
[16]:
sim.plot_variance('x')
[17]:
sim.plot_average('work')
[18]:
sim.plot_average('heat')
Plot the PDF of the position
[19]:
sim.animate_pdf('x', x_range = [-5, 5], y_range= [0, 0.5])
Plot the PDF of work
[20]:
sim.animate_pdf('work', x_range = [-2, 7], y_range = [0, 1.5])
To perform simulation swith a different potential, a new simulator has to be created. Let us create one for a non harmonic potential
[21]:
def double_well_potential(x,t):
"""Quartic potential with a double well with moving center"""
tf = 5.0
v = 0.5
if t < tf:
return (x-v*t)**4 - (x-v*t)**2
if t >= tf:
return (x-v*tf)**4 - (x-v*tf)**2
simulator2 = Simulator(
tot_sims = 100_000,
dt = 0.001,
tot_steps = 10_000,
snapshot_step = 100,
harmonic_potential= False,
potential = double_well_potential
)
[22]:
simulator2.run()
[23]:
sim2 = simulator2.simulation[0]
[24]:
sim2.animate_pdf('x', x_range=[-2, 4.5], y_range = [0, 1])